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A  number  of  magnetic  anomalies  are  present  along  the  northern  edge  of  the  lunar  South  Pole-Aitken 
(SPA)  basin.  A  variety  of  hypotheses  for  their  formation  have  been  proposed,  but  an  in-depth  study  of 
their  properties  has  not  been  performed.  Here  we  use  two  different  methods  to  invert  for  their  source 
body  characteristics:  one  that  completely  searches  a  small  parameter  space  of  less  than  ten  uniform- 
strength  dipoles  per  anomaly,  and  another  that  uses  grids  of  hundreds  of  dipoles  with  variable  magneti¬ 
zation  strengths.  Both  methods  assume  uniform  magnetization  directions  at  each  anomaly  and  with  one 
exception,  produce  nearly  the  same  results.  We  introduce  new  Monte  Carlo  methods  to  quantify  errors 
in  our  inversions  arising  from  Gaussian  time-dependent  changes  in  the  external  field  and  the  uncertain 
geometry  of  the  source  bodies.  We  find  the  errors  from  uncertainty  in  source  body  geometry  are  almost 
always  higher.  We  also  find  a  diverse  set  of  magnetization  directions  around  SPA,  which  we  combine 
with  other  physical  arguments  to  conclude  that  the  source  bodies  were  likely  magnetized  in  a  dynamo 
field.  Igneous  intrusions  are  a  reasonable  explanation  (Purucker  et  al.,  2012)  for  the  directional  variabil¬ 
ity,  since  they  could  be  intruded  over  different  magnetic  epochs.  However,  the  directional  variability  also 
implies  either  surprisingly  large  amounts  of  true  polar  wander  or  a  dynamo  not  aligned  with  the  lunar 
spin  axis.  We  also  explore  the  possibility  that  true  polar  wander  caused  by  the  SPA  impact  could  allow 
iron-rich  SPA  ejecta  to  record  a  diverse  set  of  magnetic  field  directions.  Some  of  this  material  may  have 
also  become  “sesquinary”  ejecta  and  re-impacted  across  the  Moon  on  104-106  year  timescales  to  cap¬ 
ture  such  changes.  No  completely  satisfactory  answer  emerges,  except  that  the  dipole-axis  of  the  lunar 
dynamo  may  have  been  variable  in  direction. 

Published  by  Elsevier  Inc. 


1.  Introduction 

The  South  Pole-Aitken  (SPA)  Basin  is  the  Moon’s  largest  and 
oldest  well-defined  basin  (Garrick-Bethell  and  Zuber,  2009).  It  is 
also  the  site  of  the  largest  grouping  of  magnetic  anomalies  on  the 
Moon  (Purucker,  2008).  Initially  discovered  by  the  Apollo  15  sub¬ 
satellite  (Coleman  et  al.,  1972),  the  origin  of  these  features  remains 
unknown.  Over  the  last  20  years,  three  dominant  hypotheses  for 
their  formation  have  emerged.  The  first  is  that  ejecta  from  the  Im- 
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brium  and  Serenitatis  basins  has  accumulated  at  their  antipodes, 
which  fall  within  the  SPA  basin  (Hood  and  Huang,  1991  ;Lin  et  al., 
1998;  black  crosses  in  Fig.  1).  This  antipodal  ejecta  then  attains  a 
remanent  magnetization  from  either  a  dynamo  field,  or  a  field  re¬ 
lated  to  the  impact  event  (Hood  and  Artemieva,  2008).  The  second 
hypothesis  is  that  the  anomalies  arise  from  subsurface  dikes  that 
cooled  in  a  dynamo  field  (Purucker  et  al.,  2012).  The  third  is  that 
iron-rich  material  derived  from  the  SPA  impactor  cooled  in  a  dy¬ 
namo  field  (Wieczorek  et  al.,  2012).  Determining  which,  if  any,  of 
these  hypotheses  are  true  would  have  implications  for  understand¬ 
ing  the  nature  and  history  of  the  lunar  dynamo. 

Despite  the  importance  of  these  anomalies  in  understanding  lu¬ 
nar  magnetism,  no  detailed  studies  of  the  source  body  character¬ 
istics  have  been  performed.  Purucker  et  al.  (2012)  modeled  one  of 
the  anomalies  within  SPA  as  a  series  of  vertically  magnetized  dikes, 
and  estimated  their  magnetization  strength.  Global  studies  of  lu¬ 
nar  magnetic  anomalies  have  neglected  the  SPA  region  (Arkani- 
Hamed  and  Boutin,  2014;  Takahashi  et  al.,  2014),  because  of  the 
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Fig.  1.  Lunar  South  Pole-Aitken  basin  study  areas.  The  SPA  rim  (thick  black  line)  is  from  Garrick-Bethell  and  Zuber  (2009).  (a)  Magnetic  contours  superimposed  over  LOLA 
topography  data  (contours  are  taken  from  the  map  in  part  b).  Contours  shown  are  5,  8  and  11  nT;  (b)  magnetic  field  map  with  study  areas  highlighted  (dashed  black  lines). 
Magnetic  field  data  are  from  the  Kaguya  spacecraft  magnetometer  measurements  taken  on  day  96  in  2009;  at  a  mean  altitude  of  39.4  km.  Black  crosses  indicate  the  location 
of  the  Imbrium  (west)  and  Serenitatis  (east)  basin  antipodes,  (c)  Magnetic  field  map  derived  from  a  spherical  harmonic  model  evaluated  at  39.4  km  (Purucker  and  Nicholas, 
2010).  Study  areas  are  highlighted.  We  refer  to  areas  1-3  as  the  northwestern  cluster,  and  areas  5-8  as  the  eastern  stripes. 


complicated  field  structure  in  the  region.  Using  a  method  they  re¬ 
fer  to  as  surface  vector  mapping,  Tsunakawa  et  al.  (2014)  calculated 
the  statistics  of  the  declination  of  the  field  over  the  anomalies  at 
SPA.  Based  on  the  distribution  of  declinations,  they  found  evidence 
that  the  source  bodies  are  horizontally  elongated  in  the  east-west 
direction.  While  they  also  estimated  the  magnetization  direction 
and  magnetic  paleopole  from  the  anomaly  centered  on  the  Leibnitz 
crater  within  SPA,  they  did  not  estimate  the  source  body  character¬ 
istics  of  any  of  the  larger  magnetic  anomalies  that  characterize  the 
region. 

Here,  we  will  show  that  many  individual  anomalies  can  be 
identified  within  SPA.  We  use  software  that  can  readily  generate 


and  compare  magnetic  field  maps  from  all  available  observations, 
helping  us  avoid  spurious  or  poorly  defined  magnetic  anomalies 
that  may  have  complicated  other  efforts  to  study  the  region.  Com¬ 
parisons  with  a  spherical  harmonic  model  of  the  magnetic  field 
(Purucker  and  Nicholas,  2010)  provide  a  further  test  of  consistency. 
We  will  show  that  the  source  body  magnetization  directions  offer 
a  key  test  for  their  formation  hypotheses,  and  we  make  a  substan¬ 
tial  effort  to  characterize  the  uncertainty  in  these  directions.  The 
paper  is  organized  as  follows:  Section  2  presents  our  methods,  in¬ 
cluding  uncertainty  estimation,  Section  3  presents  our  results,  and 
in  Section  4,  we  discuss  possible  origins  of  the  observed  diver¬ 
sity  in  magnetization  directions.  We  consider  impact  and  dynamo 
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origins  for  these  anomalies,  and  attempt  to  explain  their  direc¬ 
tional  variability  with  secular  variation,  a  dynamo  not  aligned  with 
the  lunar  spin  axis,  and  both  long  and  short-timescale  true  polar 
wander.  Finally,  Section  5  contains  our  conclusions. 

2.  Methods 

2.2.  Data  sources 

We  use  magnetometer  data  from  two  independent  sources:  Lu¬ 
nar  Prospector  (LP-MAG)  and  SELENE/Kaguya  (K-MAG).  The  Lu¬ 
nar  Prospector  fluxgate  magnetometer  measured  the  vector  mag¬ 
netic  field  at  up  to  18  Hz  and  transmitted  its  measurements  at 
a  reduced  resolution  of  9 Hz.  Level  IB  LP-MAG  data  are  ob¬ 
tained  from  the  NASA  Planetary  Data  System  (ppi.pds.nasa.gov). 
The  SELENE/Kaguya  spacecraft  used  a  tri-axial  fluxgate  mag¬ 
netometer  with  a  sampling  rate  up  to  32  Hz.  K-MAG  mag¬ 
netometer  data  are  obtained  from  the  SELENE  data  archives 
(12db.selene.darts.isas.jaxa.jp).  The  cadence  of  measurements  used 
in  this  study  is  0.2  Hz  for  LP-MAG  and  0.25  Hz  for  K-MAG.  Topogra¬ 
phy  data  are  from  the  Lunar  Orbiter  Laser  Altimeter  (LOLA)  (Smith 
et  al.,  2010)  (pds-geosciences.wustl.edu). 

To  best  capture  the  Moon’s  crustal  field,  all  data  used  for  anal¬ 
ysis  were  collected  in  either  the  lunar  wake  or  while  the  Moon 
was  in  the  Earth’s  magnetotail  (hereafter  referred  to  as  wake  and 
tail  data),  avoiding  distortions  caused  by  the  solar  wind  noted  by 
Halekas  et  al.  (2008);  Kurata  et  al.  (2005).  Tail  datasets  specifi¬ 
cally  exclude  epochs  during  which  plasma  sheet  disturbances  were 
noted  (Halekas  et  al.,  2012).  Consecutive  orbits  are  ~1°  in  longi¬ 
tude  apart.  At  0.2  Hz,  successive  magnetometer  measurements  are 
separated  by  ~8km  in  the  latitudinal  direction.  All  data  used  in 
this  study  are  from  the  final  months  of  Lunar  Prospector  in  1999 
and  Kaguya  in  2009,  when  measurements  were  taken  at  observa¬ 
tion  altitudes  below  50  km. 

2.2.  Data  processing  and  anomaly  identification 

Details  on  the  generation  of  magnetic  field  maps  used  in  this 
paper  (e.g.  Fig.  8)  are  after  Hemingway  and  Garrick-Bethell  (2012), 

summarized  here  for  completeness.  After  subtracting  the  back¬ 
ground  field  (taken  to  be  the  mean  field  across  each  orbit  segment 
spanning  the  region  of  interest),  the  remaining  fields  are  assumed 
to  be  due  to  crustal  sources.  Data  from  consecutive  orbits  are  com¬ 
bined  and  fit  to  square  meshes  using  Delaunay  triangulation,  with 
a  grid  spacing  of  0.25°  x  0.25°  (7.6  km  x  7.6  km  equatorial).  This 
spacing  is  finer  than  the  spacing  between  observations  (e.g.,  ~8  km 
latitude,  ~30km  longitude  for  LP-MAG),  ensuring  no  loss  of  sig¬ 
nal  variation  during  grid  generation.  Orbit  segments  that  visually 
appear  distorted  by  transient  magnetic  field  oscillations,  despite 
these  efforts,  are  discarded.  Gaps  in  spacecraft  coverage,  if  any,  can 
be  seen  in  magnetic  field  maps,  as  the  spacecraft  measurement  lo¬ 
cations  are  shown  in  all  figures  (e.g.  Fig.  8,  white  points). 

To  analyze  the  northern  SPA  region,  we  divide  the  strongest 
magnetic  anomalies  into  ten  study  areas  (Fig.  1).  We  choose 
anomalies  whose  Kaguya  and  LP  magnetic  field  maps  are  con¬ 
sistent  with  maps  from  a  spherical  harmonic  model  of  the  field 
(Purucker  and  Nicholas,  2010),  which  can  be  seen  by  comparing 
Fig.  lb  and  c.  Anomalies  that  were  not  consistent  with  this  model 
or  showed  artifacts  of  external  disturbances  were  not  used.  Areas 
1  and  2  encompass  two  approximately  linear  and  perpendicular 
magnetic  features;  in  their  vicinity  is  area  3,  also  in  proximity  to 
the  Van  de  Graaff  crater  (Dyal  et  al.,  1974).  Area  1  is  very  similar 
to  the  area  modeled  by  Purucker  et  al.  (2012).  Together,  we  refer  to 
Areas  1-3  as  the  northwestern  cluster.  Area  9  is  to  the  southwest 
of  these  areas  and  is  associated  with  the  “swirl”  albedo  anomalies 
at  Mare  Ingenii  (Blewett  et  al.,  2011;  Kramer  et  al.,  2011).  Area  4 


is  an  isolated  anomaly  to  the  southeast  of  areas  1-3.  Areas  5-8 
are  roughly  east-west  trending  linear  features  found  to  the  east  of 
these  study  areas,  and  together  we  refer  to  them  as  the  eastern 
stripes.  The  anomalies  at  areas  1  and  8  appear  connected  in  some 
maps  of  the  magnetic  field  (Purucker  et  al.  2012).  Finally,  we  also 
define  area  10,  interesting  because  it  falls  outside  the  basin’s  outer 
topographic  rim  (Garrick-Bethell  and  Zuber,  2009). 

We  divide  all  data  in  this  paper  by  the  day  of  observation,  such 
that  each  dataset  for  a  given  anomaly  is  collected  at  approximately 
the  same  altitude.  The  datasets  may  be  referenced  to  the  day 
on  which  they  were  collected,  e.g.,  1999  day  172  (LP-MAG;  read 
as  99172)  or  2009  day  123  (K-MAG;  read  as  09123).  Figs.  2  and 
3  show  an  overview  of  the  total  magnetic  field  maps  created  for  all 
observation  days  used  for  all  study  areas;  Appendix  A  shows  de¬ 
tails  of  magnetic  inversions  for  all  these  areas  (see  Sections  2.3  and 
2.4).  Supplementary  Appendix  B  shows  maps  of  the  total  field  for 
all  available  days  that  were  considered,  but  not  used. 


2.3.  Inversion  algorithm  I:  defined  dipoles ,  constant  magnetization 
(DD-CM) 

To  assess  the  robustness  of  our  results,  we  use  two  different  al¬ 
gorithms  to  invert  for  source  body  characteristics.  Both  are  regres¬ 
sions  to  find  the  least  squared  error.  The  first  completely  searches 
a  parameter  range  using  manually  placed  dipoles  as  the  source 
bodies.  The  regression  parameter  space  varies  burial  depth  in  km, 
magnetic  dipole  moment  in  Am2,  magnetic  dip  (inclination)  in  de¬ 
grees  downward  from  the  horizontal,  and  declination  in  degrees 
clockwise  from  north.  All  dipoles  in  a  given  study  area  are  con¬ 
strained  to  be  at  the  same  depth,  moment  and  direction.  Source 
dipole  positions  (latitude  0S  and  longitude  0S)  are  placed  manu¬ 
ally  by  the  user  at  locations  where  the  magnetic  field  strength  is 
greatest  according  to  our  maps.  We  refer  to  this  henceforth  as  the 
defined  dipoles,  constant  magnetization  (DD-CM)  algorithm.  As  we 
will  show,  the  ability  to  completely  search  a  parameter  space  at 
relatively  fast  computation  speeds  gives  this  algorithm  some  ad¬ 
vantages  when  uncertainties  must  be  estimated. 

The  depth  is  allowed  to  vary  between  1  km  and  99  km  in  steps 
of  0.25  km,  and  the  dipole  magnetic  moment  between  1011  Am2 
and  1014Am2  in  steps  of  2.5xl011Am2.  Depth  solutions  are  fur¬ 
ther  constrained  to  be  no  shallower  than  the  deepest  negative  to¬ 
pography  in  a  given  study  area,  ensuring  that  all  solutions  lie  be¬ 
neath  the  lunar  surface;  depth  is  measured  against  a  reference 
sphere  of  1737.4  km  (Smith  et  al.,  2010).  The  magnetic  inclination 
is  allowed  to  vary  from  -90°  to  +90°  in  steps  of  1°  and  declination 
from  -180°  to  +180°  in  steps  of  1°. 

For  n  field  measurements  inside  a  given  study  area,  the  differ¬ 
ence  between  the  model  and  the  data  is  computed.  We  minimize 
the  root  mean  square  (RMS)  total  error  of  the  Ith  east  ( SBeast ),  north 
(SB north)  and  radial  (8Bmdiai)  component  errors,  SBtotah  as: 


(SBtotai)ms  —  -  ^  (8Bi 
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- SB : 
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This  quantity  is  then  ranked  in  decreasing  order  of  total  error 
to  find  the  best-fit  magnetic  characteristics. 


2.4.  Inversion  algorithm  2:  gridded  dipoles ,  variable  magnetization 
(GD-VM) 

We  test  the  validity  of  our  inversion  results  using  an  alterna¬ 
tive,  more  objective  approach,  at  the  expense  of  exploring  a  much 
larger  parameter  space.  Instead  of  placing  dipoles  manually,  we 
establish  a  0.25°  x  0.25°  grid  of  dipoles  covering  the  entire  study 
area  (Hemingway  and  Garrick-Bethell,  2012;  Nicholas  et  al.,  2007). 
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-168  -166  -164  -162  -160 
09151,  13.2  km 


-184  -182 

99142,  21.8  km 


-18 


-192  -188  -184 
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09096,  37.6  km 


-186  -184  -182 
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Fig.  2.  Overview  of  ail  datasets  used  for  SPA  areas  1-5.  Maps  show  total  field  strength  in  units  of  nanotesla,  as  a  function  of  latitude  and  longitude  (degrees).  The  number 
in  parentheses  denotes  a  study  area  from  Fig.  1.  Mean  measurement  altitudes  are  shown.  Data  collection  days  may  be  read  in  the  format  YY-DDD,  e.g.,  09123  is  2009  day 
123,  and  so  on. 


All  dipoles  are  again  constrained  to  be  at  the  same  depth,  inclina¬ 
tion  and  declination.  Unlike  the  first  algorithm,  however,  the  mag¬ 
netic  moment  is  allowed  to  vary  among  the  dipoles.  This  creates 
a  large  m  +  3  parameter  space,  where  m  is  the  number  of  dipoles; 
Table  1  lists  the  number  of  gridded  dipoles  for  each  study  area, 
which  ranges  from  320-1408.  The  solution  is  then  found  via  a  ge¬ 
netic  search  algorithm  that  minimizes  the  RMS  error  (for  details, 
see  Hemingway  and  Garrick-Bethell,  2012).  The  algorithm  evolves 
through  “generations”  to  progress  toward  a  better  fit  (smaller  er¬ 
ror)  to  the  data.  However,  the  evolutionary  nature  of  the  algo¬ 
rithm  does  not  guarantee  optimality  of  the  solution.  We  refer  to 
this  henceforth  as  the  gridded  dipoles,  variable  magnetization  (GD- 
VM)  algorithm.  The  grid  of  best-fit  magnetic  moments  found  by 
this  method,  ranging  across  several  orders  of  magnitude  at  each 
anomaly,  can  be  seen  in  Supplementary  Appendix  C,  which  shows 
representative  results  for  each  study  area. 

2.5.  Uncertainty  estimation  for  magnetization  directions 

Error  in  our  regressions  arises  from:  ( 1 )  time-variable  contribu¬ 
tions  to  the  measured  field  due  to  non-crustal  fields  or  instrument 
noise,  (2)  the  ideal  assumption  that  the  source  is  uniformly  mag¬ 
netized  in  the  same  direction,  and  ( 3 )  the  simplified  geometry  of 
our  source  models,  even  if  the  assumption  of  unidirectional  mag¬ 
netization  were  completely  true. 

2.5.1.  Time-variable  fields 

One  method  to  account  for  the  first  source  of  error  is  to  re¬ 
port  the  set  of  solutions  with  an  RMS  error  equal  to  or  less  than 
the  uncertainty  in  the  magnetic  field  measurements  (Parker,  1991). 


The  range  of  magnetization  directions  in  this  set  then  defines  the 
error  in  the  best-fit  direction.  The  drawback  of  this  method  is  that 
the  use  of  the  field’s  uncertainty  is  arbitrary.  For  example,  at  dif¬ 
ferent  observation  altitudes,  the  magnitude  of  the  RMS  error  for 
the  best-fit  solution  will  vary,  such  that  the  measurement  uncer¬ 
tainty  will  have  a  different  impact  on  the  size  of  the  set  of  solu¬ 
tions  that  represent  the  error.  More  importantly,  as  we  show  be¬ 
low,  the  set  of  allowable  solutions  suggested  by  propagating  Gaus¬ 
sian  measurement  noise  (or  external  field  oscillations)  of  a  given 
magnitude  through  the  inversion  is  often  small  compared  with  the 
size  of  the  set  of  solutions  contained  by  an  RMS  error  of  equal 
magnitude.  That  is,  using  the  uncertainty  in  the  field  measurement 
to  define  the  allowable  solution  set  can  overestimate  the  error  in 
the  regression. 

To  demonstrate  this  effect,  we  simulate  the  effect  of  Gaus¬ 
sian  measurement  noise  in  each  of  the  45  regressions  we  perform 
(all  maps  shown  in  Figs.  2  and  3).  Before  doing  so,  we  first  es¬ 
timate  the  characteristic  noise  magnitude.  For  both  LP-MAG  and 
K-MAG,  the  instrument  noise  is  <0.1  nT,  such  that  the  dominant 
source  of  error  arises  from  fluctuations  in  the  interplanetary  mag¬ 
netic  field  (IMF)  rather  than  the  instrument.  To  estimate  a  charac¬ 
teristic  amplitude  of  the  oscillations  in  the  IMF  over  the  timescale 
of  observation  at  one  magnetic  anomaly,  while  in  the  wake  or 
magnetotail,  one  would  ideally  examine  the  time-dependent  oscil¬ 
lations  just  prior  to  flying  over  the  anomaly  of  interest.  However, 
there  is  no  simple  way  to  separate  oscillations  arising  from  small- 
scale  crustal  fields  and  time  dependent  fields.  As  a  rough  proxy, 
we  can  examine  the  field  variations  over  the  Mare  Imbrium  region 
while  the  Moon  is  in  the  wake.  Since  the  Imbrium  region  is  known 
to  have  the  lowest  amount  of  crustal  magnetism  anywhere  on  the 
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Fig.  3.  Overview  of  all  datasets  used  for  SPA  areas  6-10.  Maps  show  total  field  strength  in  units  of  nanotesla,  as  a  function  of  latitude  and  longitude  (degrees).  The  number 
in  parentheses  denotes  a  study  area  from  Fig.  1.  Mean  measurement  altitudes  are  shown.  Data  collection  days  may  be  read  in  the  format  YY-DDD,  e.g.,  09123  is  2009  day 
123,  and  so  on. 


Moon  (Mitchell  et  al.,  2008;  Tsunakawa  et  al.,  2014),  magnetic  field 
oscillations  there  are  most  likely  to  represent  only  time  dependent 
effects.  We  find  that  the  standard  deviations  of  the  field  at  Im- 
brium  are  0.07  nT  and  0.13  nT  for  two  representative  days  (Fig.  4). 
However,  to  be  conservative,  we  use  an  order  of  magnitude  higher 
standard  deviation  of  1  nT  as  the  characteristic  value.  Takahashi  et 
al.  (2014)  make  a  similar  (1  nT)  assumption  about  the  magnitude 
of  external  field  contributions. 

Next,  to  simulate  the  effect  of  typical  time-dependent  fields  in 
our  DD-CM  regressions,  we  added  random  Gaussian  noise  with  a 
standard  deviation  of  1  nT  to  each  field  component  (east,  north, 
radial)  for  all  45  datasets.  An  example  is  shown  in  Fig.  5.  Note  that 
we  are  assuming  that  the  noise  at  each  measurement  point  is  un¬ 
correlated  with  the  previous  measurement  point.  Our  assumption 
is  only  valid  if  oscillations  in  the  IMF  occur  on  timescales  faster 
than  the  time  between  measurements  (5  seconds  for  LP  data,  and 
4  seconds  for  Kaguya  data).  Future  work  will  determine  the  depen¬ 
dence  of  our  error  estimates  on  a  variety  of  correlation  timescales, 
and  better  determine  the  typical  timescale  of  time-dependent  field 
oscillations. 

For  each  of  the  45  datasets  (Figs.  2  and  3),  we  generated  20 
sets  of  noise-added  data.  Assuming  the  best-fit  directions  follow 
a  Fisher  distribution  (the  Fisher  distribution  gives  the  probability 
of  finding  a  direction  within  an  angular  area  from  the  true  mean; 


Fisher,  1953),  their  angular  standard  deviations  and  Fisher  distri¬ 
bution  precision  parameters  /<,  were  then  obtained  from  each  set 
of  20  simulations.  The  mean  k  value  across  all  altitudes  within  a 
study  area  was  then  used  as  the  final  estimate  of  the  uncertainty 
from  time-dependent  fields  for  that  study  area  (for  example,  for 
area  1,  which  uses  data  from  four  different  altitudes,  four  differ¬ 
ent  k  values  were  averaged  together).  The  ability  to  simulate  the 
effect  of  noise  in  this  manner  represents  an  advantage  of  the  DD- 
CM  method,  which  has  a  faster  computation  time  compared  to  the 
GD-VM  method. 

We  find  that  the  angular  standard  deviation  in  best-fit  direc¬ 
tion  associated  with  1  nT  noise  is  always  smaller  than  the  set  of 
directions  permitted  when  using  an  error  threshold  of  1  nT.  For 
example,  the  best-fit  RMS  error  for  area  1,  09096  (Fig.  5)  is  2.4  nT, 
which  is  already  greater  than  the  ~1  nT  uncertainty  we  estimate 
for  the  field  (this  inequality  is  true  for  nearly  all  inversions  in  this 
study).  The  set  of  solutions  that  are  allowed  by  considering  1  nT 
additional  RMS  error  (total  residual  error  3.4  nT)  produces  ~30° 
of  uncertainty  in  the  best-fit  magnetization  direction  (Fig.  6).  In 
contrast,  Fig.  6  also  shows  that  adding  1  nT  Gaussian  noise  to  the 
observations  via  the  above  method  only  produces  an  angular  stan¬ 
dard  deviation  of  5°  in  the  best-fit  direction.  The  smaller  size  of 
the  effect  of  simulated  Gaussian  noise  is  true  in  all  of  our  regres¬ 
sions. 
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Table  1 

Lunar  Prospector  (LP)  and  Kaguya  (KG)  datasets  used  for  inversions,  with  best-fit  parameters  found  using  both  the  DD-CM  and  GD-VM  algorithms.  Mean  observational 
altitude,  data  day  and  site  latitude  and  longitude  are  shown.  GD-VM  dipoles  column  indicates  the  number  of  gridded  dipoles  used  for  the  GD-VM  simulation.  Estimated 
magnetizations  from  both  algorithms  are  shown  (see  text  for  details).  The  merged  row  for  each  study  area  shows  the  best-fit  parameters  obtained  from  merged  data  at  all 
altitudes  using  the  GD-VM  algorithm.  Highlighted  rows  correspond  to  figures  shown  in  the  main  text  (Figs.  8-17);  GD-VM  dipole  maps  for  these  datasets  are  also  shown  in 
Supplementary  Appendix  C.  Best-fit  figures  for  non-highlighted  rows  are  in  Appendix  A. 


DD-CM  algorithm  GD-VM  Algorithm 


Study  Area 

Spacecraft 

Year 

Day 

Latitude 
(deg  N) 

Longitude 
(deg  E) 

Ait 

(km) 

Depth 

(km) 

Moment 

(1012  Am2) 

Inc 

(deg) 

Dec 

(deg) 

Mag 

(A/m) 

Dipoles 

Depth 

(km) 

Inc 

(deg) 

Dec 

(deg) 

Mag 

(A/m) 

1 

LP-MAG 

1999 

142 

-25 

-184 

22.4 

25 

8.4 

13 

-123 

0.12 

580 

5 

18 

-128 

0.70 

LP-MAG 

1999 

61 

-25 

-184 

32.3 

24 

6.4 

32 

-116 

0.09 

580 

5 

22 

-128 

0.62 

K-MAG 

2009 

96 

-25 

-184 

39.4 

36 

11.0 

42 

-132 

0.06 

580 

26 

30 

-126 

0.12 

K-MAG 

2009 

123 

-25 

-184 

33 

37 

6.4 

20 

-114 

0.11 

580 

14 

20 

-116 

0.24 

Merged 

580 

6 

19 

-123 

0.42 

2 

LP-MAG 

1999 

33 

-20 

-188 

32.8 

46 

21.0 

36 

21 

0.06 

1408 

21 

11 

-71 

0.20 

K-MAG 

2009 

69 

-20 

-188 

44.7 

37.5 

26.0 

29 

-129 

0.06 

1408 

39 

5 

-54 

0.10 

K-MAG 

2009 

96 

-20 

-188 

39.4 

40 

15.0 

26 

-135 

0.07 

1408 

38 

7 

-52 

0.11 

K-MAG 

2009 

123 

-20 

-188 

33 

40 

20.0 

32 

-141 

0.07 

1408 

26 

7 

-53 

0.16 

Merged 

1408 

12 

8 

-59 

0.14 

3 

LP-MAG 

1999 

61 

-27 

-192 

32.3 

47 

19.0 

31 

1 

0.08 

432 

21 

20 

3 

0.16 

K-MAG 

2009 

69 

-27 

-192 

44.7 

55 

62.0 

32 

-8 

0.10 

432 

27 

13 

-1 

0.16 

K-MAG 

2009 

96 

-27 

-192 

39.4 

49.5 

18.0 

49 

4 

0.07 

432 

23 

27 

-8 

0.13 

K-MAG 

2009 

123 

-27 

-192 

33 

48 

17.0 

49 

2 

0.06 

432 

23 

21 

2 

0.13 

Merged 

432 

11 

25 

-8 

0.15 

4 

LP-MAG 

1999 

115 

-33 

-183 

32.3 

28 

6.4 

-64 

-46 

0.06 

320 

24 

-55 

-56 

0.09 

LP-MAG 

1999 

142 

-33 

-183 

22.4 

40.5 

12.0 

-73 

-89 

0.07 

320 

13 

-52 

-50 

0.20 

K-MAG 

2009 

96 

-33 

-183 

39.4 

28 

6.2 

-69 

-23 

0.06 

320 

24 

-70 

-28 

0.08 

K-MAG 

2009 

151 

-33 

-183 

13.7 

23.5 

5.7 

-82 

86 

0.06 

320 

8 

-62 

-42 

0.30 

Merged 

320 

9 

-85 

-73 

0.17 

5 

K-MAG 

2009 

96 

-19 

-165 

39.4 

32.5 

11.0 

-33 

-16 

0.06 

850 

24 

-23 

-23 

0.15 

K-MAG 

2009 

151 

-19 

-165 

13.7 

55 

34.0 

-24 

-12 

0.12 

850 

25 

-34 

-33 

0.17 

LP-MAG 

1999 

142 

-19 

-165 

22.4 

46 

21.0 

-29 

-23 

0.09 

850 

26 

-30 

-33 

0.14 

K-MAG 

2009 

123 

-19 

-165 

33 

38.5 

7.9 

-56 

6 

0.06 

850 

28 

-36 

-25 

0.11 

Merged 

850 

12 

-30 

-39 

0.14 

6 

LP-MAG 

1999 

142 

-24 

-165 

22.4 

46 

23.0 

47 

-38 

0.10 

576 

23 

50 

-28 

0.16 

LP-MAG 

1999 

33 

-24 

-165 

32.8 

55 

33.0 

49 

24 

0.08 

576 

31 

39 

53 

0.13 

LP-MAG 

1999 

115 

-24 

-165 

32.3 

55 

30.0 

41 

-8 

0.08 

576 

21 

47 

18 

0.16 

K-MAG 

2009 

96 

-24 

-165 

39.4 

46 

22.0 

49 

31 

0.07 

576 

32 

31 

-21 

0.12 

K-MAG 

2009 

123 

-24 

-165 

33 

55 

21.0 

43 

-14 

0.08 

576 

34 

38 

-25 

0.11 

Merged 

576 

12 

43 

15 

0.14 

7 

LP-MAG 

1999 

33 

-27 

-165 

32.8 

55 

22.0 

35 

-160 

0.09 

341 

26 

28 

-132 

0.15 

LP-MAG 

1999 

115 

-27 

-165 

32.3 

41.5 

12.0 

45 

-180 

0.07 

341 

25 

38 

159 

0.12 

K-MAG 

2009 

96 

-27 

-165 

39.4 

24.5 

8.0 

46 

152 

0.08 

341 

21 

34 

155 

0.14 

K-MAG 

2009 

123 

-27 

-165 

33 

37 

11.0 

40 

131 

0.07 

341 

24 

34 

143 

0.12 

Merged 

341 

14 

39 

165 

0.13 

8 

K-MAG 

2009 

123 

-31 

-164 

33 

25.5 

6.9 

2 

146 

0.05 

1100 

39 

-3 

-131 

0.08 

LP-MAG 

1999 

61 

-31 

-164 

32.3 

46 

21.0 

1 

-172 

0.03 

1100 

45 

3 

134 

0.06 

K-MAG 

2009 

151 

-31 

-164 

13.7 

25 

5.9 

-17 

138 

0.05 

1100 

13 

0 

123 

0.31 

LP-MAG 

1999 

142 

-31 

-164 

22.4 

55 

20.0 

-17 

-128 

0.04 

1100 

36 

-6 

-119 

0.10 

K-MAG 

2009 

96 

-31 

-164 

39.4 

28 

7.8 

-1 

-162 

0.03 

1100 

46 

-3 

-132 

0.07 

Merged 

1100 

13 

-1 

-119 

0.12 

9 

LP-MAG 

1999 

61 

-36 

-198 

32.3 

29.5 

26.0 

1 

119 

0.16 

672 

7 

-4 

145 

0.52 

LP-MAG 

1999 

88 

-36 

-198 

32.2 

34 

37.0 

-12 

146 

0.13 

672 

9 

-11 

162 

0.43 

LP-MAG 

1999 

115 

-36 

-198 

32.3 

25 

30.0 

-6 

139 

0.14 

672 

7 

-8 

141 

0.49 

LP-MAG 

1999 

142 

-36 

-198 

22.4 

22.5 

17.0 

-2 

144 

0.18 

672 

4 

-5 

153 

0.95 

K-MAG 

2009 

96 

-36 

-198 

39.4 

26 

23.0 

-5 

125 

0.16 

672 

7 

-8 

141 

0.41 

K-MAG 

2009 

123 

-36 

-198 

33 

23.5 

15.0 

-5 

135 

0.19 

672 

5 

-3 

143 

0.65 

Merged 

672 

4 

-7 

146 

0.58 

10 

LP-MAG 

1999 

61 

-14 

-153 

32.3 

36 

25.0 

34 

-180 

0.05 

506 

15 

15 

156 

0.20 

LP-MAG 

1999 

88 

-14 

-153 

32.2 

20.5 

8.3 

32 

166 

0.10 

506 

12 

11 

-167 

0.23 

LP-MAG 

1999 

142 

-14 

-153 

22.4 

21 

17.0 

29 

-175 

0.06 

506 

5 

14 

-174 

0.56 

K-MAG 

2009 

123 

-14 

-153 

33 

19 

15.0 

21 

-172 

0.06 

506 

5 

10 

-172 

0.50 

K-MAG 

2009 

151 

-14 

-153 

13.7 

16 

14.0 

25 

-168 

0.07 

506 

4 

10 

-155 

0.87 

Merged 

506 

4.5 

13 

-160 

0.47 

2.5.2.  Non-uniform  magnetization  direction 

The  second  contribution  to  error  is  the  non-uniform  magneti¬ 
zation  directions  in  the  source  bodies.  However,  it  is  impossible 
to  make  any  inferences  about  the  magnetization  if  we  permit  the 
infinite  number  of  possible  source  magnetizations  with  mixed  di¬ 
rections.  Therefore  we  must  at  least  assume  that  the  source  is  uni¬ 
formly  magnetized.  We  attempt  to  mitigate  this  effect  by  selecting 
small,  well-defined  anomalies  (Fig.  1).  Of  course,  if  non-uniformity 
dominates  the  source  of  error  then  these  regression  results  are  less 
meaningful. 


2.5.3.  Source  body  geometry 

Finally,  we  address  what  is  likely  the  dominant  source  of  uncer¬ 
tainty:  the  complex  geometry  of  the  magnetic  source  bodies  com¬ 
pared  to  our  simplified  models.  We  do  this  in  two  ways.  The  first 
is  to  compute  the  best-fit  magnetization  separately  for  different 
altitude  datasets  at  each  anomaly,  instead  of  fitting  to  a  merged 
altitude  dataset.  Data  from  different  altitudes  display  a  diversity 
of  magnetic  field  morphologies,  as  field  strength  decays  with  al¬ 
titude  (which  ranges  in  this  study  from  13-46  km).  Hence,  each 
altitude  will  lead  to  different  choices  of  dipole  locations  for  the 
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Longitude  (degrees  E) 


Magnetic  field  strength  (nT) 


Fig.  4.  Magnetic  field  maps  and  histograms  of  field  strength  at  eastern  Mare  Imbrium  in  the  lunar  wake.  Panels  (a)  and  (c)  are  magnetic  field  maps  at  Imbrium  using 
LP-MAG  observation  days  118  and  131  in  1999,  respectively.  White  dots  indicate  observation  locations.  The  mean  observation  altitude  for  day  118  is  31.3  km  and  29.8  km  for 
day  131.  Panels  (b)  and  (d)  are  corresponding  histograms  of  the  strength  of  the  total  magnetic  field. 
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Fig.  5.  Error  estimation  method  for  time-variable  (non-crustal)  magnetic  fields.  Kaguya  magnetic  field  data  for  area  1,  2009  day  96  (left  column),  the  effect  of  random 
noise  with  a  standard  deviation  of  1  nT  added  to  that  data  (middle  column)  and  one  of  our  Monte  Carlo  trial  regressions  for  the  depth,  moment  and  direction  values  for 
data  +  noise  with  the  DD-CM  algorithm  (right  column).  Dipole  source  geometry  is  identical  to  the  placement  in  Fig.  8  (black  dots). 


DD-CM  algorithm;  for  example,  compare  Fig.  8  to  Figures  A-l  to 
A-3  (Appendix  A).  This  variability  is  an  advantage,  as  it  allows  us 
to  probe  the  sensitivity  of  our  results  to  source  geometry,  with 
the  spread  in  the  best-fit  magnetization  directions  representing  the 
measure  of  uncertainty  (this  advantage  will  be  particularly  obvious 
for  area  2,  which  exhibits  bimodal  error  minima  in  its  error  space). 
This  is  analogous  to  the  practice  in  paleomagnetism  of  sampling 
a  single  rock  formation  multiple  times  at  different  sites  (Irving, 
1964).  Because  each  measured  direction  at  each  altitude  is  inde¬ 
pendently  calculated  we  can  also  assign  a  95%  confidence  interval 
to  the  mean  direction  from  all  altitudes  (Butler,  1998).  In  practice, 


we  combine  the  directions  obtained  from  both  the  DD-CM  and  GD- 
VM  algorithms  to  calculate  the  final  mean  and  confidence  interval, 
even  though  the  direction  measurements  for  a  given  day  are  not 
strictly  independent  across  the  two  algorithms  (see  Section  3).  In 
calculating  the  final  mean,  we  also  include  one  best-fit  direction 
obtained  from  the  merged  data  at  all  altitudes,  using  the  GD-VM 
algorithm. 

The  second  way  we  account  for  uncertainty  due  to  source  body 
geometry  is  to  explore  the  range  in  directions  that  are  returned 
by  displacing  the  nominal  DD-CM  algorithm  dipoles  for  all  45 
datasets.  We  calculate  the  effect  of  placing  dipoles  anywhere  ran- 
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Fig.  6.  RMS  error  map  for  the  magnetization  direction  for  area  1,  2009  day  96,  centered  at  0°  inclination  and  0°  declination  (whole-sphere  Mollweide  projection,  in  which 
the  southern  hemisphere  is  positive  inclination,  same  as  in  Fig.  18).  The  small  white  circular  contour  indicates  one  angular  standard  deviation  of  dispersion,  from  Monte 
Carlo  simulations  of  the  addition  of  1  nT  Gaussian  noise  to  observations  (Fig.  5).  The  larger,  outer  white  dashed  contour  indicates  the  minimum  error  solution  plus  1  nT. 
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Fig.  7.  Error  estimation  method  for  uncertainty  in  source  body  geometry.  Black 
points  are  the  nominal  dipole  locations  for  a  DD-CM  regression  at  area  1.  Circles 
represent  the  allowable  locations  of  dipoles  placed  randomly  in  our  Monte  Carlo 
simulations  and  are  0.5°  in  radius.  K-MAG  total  field  observations  in  the  lunar  wake 
on  2009  day  96.  Compare  to  Fig.  8. 

domly  on  a  0.5°-radius  circle  from  their  nominal  location  (Fig.  7). 
The  value  of  0.5°  represents  the  approximate  error  in  longitudinal 
uncertainty  in  the  anomaly  peak  field  location  (orbits  are  spaced 
by  approximately  1°  for  a  given  constant  altitude  data  set).  For 
each  data  set,  we  use  100  random  placements  of  dipoles  at  any 
location  on  the  defined  circles  and  calculate  the  100  best-fit  mag¬ 
netization  directions  and  /<  values.  The  k  values  are  then  aver¬ 
aged  across  altitudes  within  an  area  to  obtain  the  uncertainty  from 
source  geometry  for  that  area.  Again,  this  represents  an  advantage 
of  the  DD-CM  method  compared  to  the  GD-VM  method,  due  to 
the  ability  to  modify  model  dipole  placements  and  a  faster  com¬ 
putation  time.  We  will  see  this  advantage  manifest  when  explor¬ 
ing  ambiguous  regression  results  for  area  2.  Finally,  we  note  that 
for  area  4  we  use  0.25°-radius  circles,  due  to  the  close  proximity 
of  the  model  dipoles  there. 

The  dispersions  in  best-fit  directions  obtained  from  our  Monte 
Carlo  simulations  for  time-dependent  fields  (Section  2.5.1)  and 
0.5°-radius  dipole  displacements  are  then  combined  using  their  k 
values,  averaged  together  for  each  area  (Irving,  1964): 

1/^tota/  =  1  /^time-dependent,  mean  T-  1  /^source-geometry,  mean  (2) 

An  angular  standard  deviation  063  is  then  estimated  by  (Irving, 
1964): 

063^81  /k03  (3) 


Finally,  it  is  linearly  added  to  the  angular  95%  confidence  inter¬ 
val  obtained  from  the  mean  of  the  variable  altitude  results  within 
a  given  area  (conservatively  assuming  the  Monte  Carlo  error  does 
not  add  in  quadrature,  since  it  may  be  correlated). 

The  following  is  a  summary  of  our  uncertainty  estimation 
methods: 

1.  Account  for  time-dependent  external  field  contributions  and 
instrument  noise:  Perform  Monte  Carlo  simulations  for  every 
anomaly  and  altitude  with  the  addition  of  1  nT  noise,  and  ob¬ 
tain  mean  k  for  the  best-fit  directions  within  an  anomaly  area, 
using  the  DD-CM  algorithm  (Fig.  5,  Table  2). 

2.  Account  for  uncertainty  in  source  geometry:  Perform  Monte 
Carlo  simulations  for  every  anomaly  and  altitude  while  altering 
the  dipole  placement,  and  obtain  mean  k  for  the  best-fit  direc¬ 
tions  within  an  anomaly  area,  using  the  DD-CM  algorithm  (Fig. 
7,  Table  2). 

3.  Combine  the  k  values  from  steps  1  and  2  and  obtain  the  angu¬ 
lar  standard  deviation  via  Eqs.  2  and  3. 

4.  Account  for  uncertainty  in  source  geometry  by  using  variable 
altitude  data,  and  hence  variable  dipole  placement.  Obtain  a 
95%  confidence  interval  on  the  mean  direction  from  the  com¬ 
bined  best-fit  directions  from  the  DD-CM  and  GD-VM  algo¬ 
rithms  (Table  3). 

5.  Linearly  add  the  angular  dispersion  from  step  3  to  the  angular 
95%  confidence  interval  from  step  4  to  obtain  the  total  angular 
uncertainty  estimate  in  the  best-fit  direction. 

3.  Results 

3.1  Regression  results 

Minimum  error  magnetization  directions  are  found  for  each  Lu¬ 
nar  Prospector  and  Kaguya  dataset,  using  both  the  DD-CM  (Section 
2.3)  and  the  GD-VM  (Section  2.4)  algorithms.  Table  1  compiles  the 
values  of  these  best-fit  directions.  The  best-fit  models  are  shown 
in  Figs.  8  -17  for  study  areas  1-10.  The  figures  show  10  represen¬ 
tative  examples,  one  for  each  study  area;  Appendix  A  shows  all 
results  for  every  dataset  listed  in  Table  1.  Angular  standard  devi¬ 
ations  and  k  for  the  noise-added  and  dipole  displacement  simula¬ 
tions  are  shown  in  Table  2  for  all  datasets.  RMS  error  maps  for 
all  datasets  (similar  to  Fig.  6)  are  compiled  in  Figs.  18  and  19; 


M.  Nayak  et  al./ Icarus  286  (2017)  153-192 


161 


Table  2 

Uncertainty  estimates  for  regression  results  from  Monte  Carlo  simulations.  Minimum  RMS  error  recovered  using  the  DD-CM  and  GD-VM 
algorithms  for  all  areas  and  altitudes.  Fisher  precision  parameter  k,  angular  standard  deviation  063  (^  81//<1/2),  obtained  from  Monte  Carlo 
simulations  for  the  effects  of  time-dependent  fields  (“noise”)  and  0.5°-radius  displaced  dipoles,  for  the  DD-CM  algorithm.  Highlighted  rows 
correspond  to  figures  shown  in  the  main  text  (Fig.  8-17). 


Study  Area 

Spacecraft 

Year 

Day 

Min.  RMS  error  Min.  RMS  error  Noise  simulation  Displaced  dipole  simulation 

(DD-CM)  (GD-VM)  063  Precision  063  Precision 

(nT)  (nT)  parameter  k  parameter  k 

1 

LP-MAG 

1999 

142 

7.4 

3.9 

1.6 

2526 

7.4 

120 

LP-MAG 

1999 

61 

5.6 

2.2 

2.1 

1474 

7.3 

124 

K-MAG 

2009 

96 

2.4 

1.1 

5.0 

264 

4.2 

381 

K-MAG 

2009 

123 

6.3 

2.3 

2.9 

774 

3.8 

446 

2 

LP-MAG 

1999 

33 

5.9 

3.1 

1.5 

3120 

4.5 

332 

K-MAG 

2009 

69 

3.7 

2.3 

3.3 

620 

5.7 

203 

K-MAG 

2009 

96 

4.0 

2.2 

2.1 

1542 

44.1 

3 

K-MAG 

2009 

123 

5.7 

3.2 

1.7 

2318 

14.3 

32 

3 

LP-MAG 

1999 

61 

4.1 

2.0 

3.1 

676 

5.7 

202 

K-MAG 

2009 

69 

3.1 

1.6 

3.1 

674 

7.2 

128 

K-MAG 

2009 

96 

1.7 

0.9 

4.3 

353 

6.1 

177 

K-MAG 

2009 

123 

2.8 

1.3 

5.3 

234 

6.4 

160 

4 

LP-MAG 

1999 

115 

1.5 

1.1 

6.4 

159 

4.7 

300 

LP-MAG 

1999 

142 

4.4 

2.0 

4.2 

375 

5.3 

232 

K-MAG 

2009 

96 

1.1 

0.8 

7.2 

127 

5.5 

219 

K-MAG 

2009 

151 

8.1 

3.7 

1.9 

1761 

8.6 

89 

5 

K-MAG 

2009 

96 

2.0 

1.0 

5.3 

232 

3.7 

472 

K-MAG 

2009 

151 

9.1 

4.0 

1.3 

3959 

3.6 

501 

LP-MAG 

1999 

142 

4.6 

1.9 

2.0 

1648 

3.3 

601 

K-MAG 

2009 

123 

2.2 

1.1 

5.2 

243 

4.1 

383 

6 

LP-MAG 

1999 

142 

5.7 

3.1 

2.4 

1150 

15.1 

29 

LP-MAG 

1999 

33 

3.8 

1.6 

4.8 

288 

11.5 

50 

LP-MAG 

1999 

115 

3.3 

2.0 

4.8 

283 

5.3 

234 

K-MAG 

2009 

96 

2.2 

1.1 

7.8 

109 

6.9 

140 

K-MAG 

2009 

123 

2.7 

1.5 

4.6 

308 

4.4 

343 

7 

LP-MAG 

1999 

33 

3.1 

1.6 

5.0 

262 

4.6 

311 

LP-MAG 

1999 

115 

2.2 

1.3 

6.0 

184 

6.3 

164 

K-MAG 

2009 

96 

1.5 

0.9 

8.9 

83 

8.3 

95 

K-MAG 

2009 

123 

2.4 

1.4 

5.2 

244 

6.7 

147 

8 

K-MAG 

2009 

123 

2.6 

1.7 

4.4 

347 

28.0 

8 

LP-MAG 

1999 

61 

2.5 

1.5 

6.0 

184 

3.9 

443 

K-MAG 

2009 

151 

7.5 

3.9 

2.4 

1121 

36.6 

5 

LP-MAG 

1999 

142 

4.4 

2.6 

3.7 

476 

3.2 

625 

K-MAG 

2009 

96 

1.6 

1.1 

7.2 

128 

5.3 

235 

9 

LP-MAG 

1999 

61 

6.3 

3.2 

1.4 

3496 

4.3 

352 

LP-MAG 

1999 

88 

5.6 

2.3 

1.8 

1999 

8.1 

99 

LP-MAG 

1999 

115 

5.7 

2.1 

1.9 

1837 

9.6 

71 

LP-MAG 

1999 

142 

9.0 

3.5 

1.0 

6997 

8.5 

91 

K-MAG 

2009 

96 

4.0 

1.5 

2.8 

853 

8.1 

101 

K-MAG 

2009 

123 

4.5 

1.5 

1.8 

2004 

5.4 

228 

10 

LP-MAG 

1999 

61 

3.0 

1.5 

3.0 

730 

7.7 

112 

LP-MAG 

1999 

88 

3.2 

1.7 

7.9 

104 

11.9 

46 

LP-MAG 

1999 

142 

4.7 

2.3 

2.2 

1315 

10.2 

63 

K-MAG 

2009 

123 

1.9 

1.3 

3.9 

431 

7.3 

124 

K-MAG 

2009 

151 

8.0 

4.2 

1.2 

4571 

15.5 

27 

Table  3 

Uncertainty  estimates  derived  from  using  different  altitude  datasets.  The  paleolatitude,  paleolon- 
gitude,  value  of  k,  and  the  95%  confidence  circle  (oc95)  are  obtained  from  the  combined  set  of 
best-fit  directions  from  the  DD-CM  and  GD-VM  algorithms,  across  all  altitudes.  The  oc95  value  is 
added  to  the  angular  uncertainty  estimated  by  Monte  Carlo  simulations  (equations  2  and  3,  de¬ 
rived  from  mean  Table  2  values)  to  obtain  the  final  circular  error,  also  shown.  The  distorted  ver¬ 
sion  of  this  final  error  circle  at  the  paleopole  is  also  shown  (Figs.  21  and  22):  the  semi-axis  of  the 
ellipse  along  the  great  circle  path  from  site  to  pole  (dp)  and  semi-axis  of  the  ellipse  perpendicular 
to  that  path  (dm). 


Area 

Altitude  averaged 

k  a -95 

Paleolatitude 

Paleolongitude 

Final 

circular  error 

dp 

dm 

1 

58.5 

6.8 

-34.8 

84.2 

12.2 

7.0 

13.1 

2 

3.3 

33.9 

8.5 

100.4 

41.0 

24.0 

44.3 

3 

38.7 

8.4 

46.9 

165.1 

15.7 

9.6 

17.4 

4 

27.8 

9.9 

51.8 

39.1 

16.4 

24.5 

28.4 

5 

33.0 

9.1 

68.1 

104.5 

13.3 

8.6 

15.1 

6 

14.5 

12.4 

38.0 

-165.1 

19.9 

16.6 

25.7 

7 

11.2 

16.1 

-78.5 

-62.0 

24.5 

18.2 

29.9 

8 

3.2 

30.9 

-54.2 

37.8 

37.1 

18.7 

37.3 

9 

50.0 

5.9 

-37.0 

-67.6 

12.6 

6.3 

12.6 

10 

25.4 

9.2 

-83.3 

82.9 

18.9 

10.3 

19.8 
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Fig.  8.  Area  1  best-fit  results  for  both  inversion  algorithms.  (Left  column)  K-MAG  observations  in  the  lunar  wake  in  2009,  day  96,  at  38.5  km,  compared  to  (middle  column)  a 
model  using  several  manually  placed  dipoles  of  equal  magnetization  (DD-CM),  and  (right  column)  a  grid-based  model  using  many  dipoles  of  variable  magnetization  (GD-VM). 
Source  dipoles  in  the  DD-CM  algorithm  are  approximately  located  at  maxima  in  the  observed  total  field  (black  dots).  White  dots  in  the  east  panel  are  locations  of  spacecraft 
magnetometer  measurements. 
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Fig.  9.  Area  2  best-fit  results  with  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  that  magnetometer  measurements  are  taken  from  K-MAG,  2009  day  123, 
altitude  32.5  km. 


these  illustrate  the  difference  between  using  an  arbitrary  uncer¬ 
tainty  threshold  and  the  Monte  Carlo  methods  employed  in  this 
work. 

Fig.  20  compiles  the  values  of  the  best-fit  directions  from  the 
DD-CM  and  GD-VM  algorithms.  For  the  DD-CM  algorithm,  we  show 
one  standard  deviation  of  dispersion  from  1  nT  noise  (Fig.  20a)  and 
one  standard  deviation  of  dispersion  from  displacing  the  model 
dipoles  by  0.5°  (Fig.  20b)  (see  Methods).  Overall,  we  find  the  error 


from  the  uncertainty  in  source  geometry  (Fig.  20b)  is  usually  larger 
than  the  error  from  the  effects  of  time-dependent  fields  (Fig.  20a). 

Interestingly,  we  find  source  geometry  errors  at  areas  2  and 
8  are  substantially  larger  than  for  the  other  anomalies.  The  large 
error  in  area  2  arises  from  the  existence  of  two  nearly  equal 
local  minima  in  its  DD-CM  error  map,  which  is  not  seen  at  any 
other  anomaly  (see  Fig.  18).  The  two  minima  are  approximately 
110°  apart.  Small  displacements  of  the  dipoles  from  their  nom- 
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Fig.  10.  Area  3  best-fit  results  with  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  the  altitude  here  is  37.6  km. 
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Fig.  11.  Area  4  best-fit  results  with  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  the  altitude  here  is  38.3  km. 


inal  locations  flip  the  best-fit  solution  into  the  other  minimum, 
producing  its  large  uncertainty  ellipses,  particularly  for  day  09096 
(Fig.  20b).  One  of  our  four  best-fit  solutions  (day  99033)  falls 
in  this  secondary  minimum.  Because  area  1  and  2  have  slight 
overlap  (Fig.  1),  we  tested  if  the  double  minima  could  be  a  result 
of  area  1  data  influencing  the  regression  at  area  2.  We  ran  the 
DD-CM  algorithm  against  data  for  area  2,  but  this  time  excluded 


spacecraft  measurements  over  area  1.  The  results  are  not  signif¬ 
icantly  different;  the  two  minima  are  still  seen  for  area  2.  Using 
the  GD-VM  algorithm,  one  finds  the  best-fit  directions  for  area  2 
are  approximately  between  the  two  minima  found  by  the  DD-CM 
algorithm  (Fig.  20c),  ~60°  from  either  one.  This  difference  between 
the  two  algorithms  is  the  largest  for  any  area.  Further  work  will 
be  required  to  better  understand  the  nature  of  the  error  space  at 
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Fig.  12.  Area  5  best-fit  results  with  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  that  magnetometer  measurements  are  taken  from  K-MAG,  2009  day  123, 
altitude  33.3  km. 
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Fig.  13.  Area  6  best-fit  results  from  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  that  magnetometer  measurements  are  taken  from  K-MAG,  2009  day  123, 
altitude  32.4  km. 


area  2,  and  its  source  body  characteristics;  we  will  investigate,  for 
example,  oppositely  magnetized  blocks  of  magnetization  in  this 
region  (cf.  Parker,  1988). 

For  area  8,  the  large  error  arises  from  the  flatness  of  the  error 
space,  instead  of  the  existence  of  multiple  minima  (Fig.  18).  In 
particular,  the  declination  of  this  nearly  equator-pointing  magneti¬ 
zation  vector  is  poorly  constrained,  with  results  spanning  a  range 
of  ~100°  of  arc  (a  result  found  to  be  true  using  either  the  DD-CM 
or  GD-VM  algorithm).  Again,  without  the  error  maps  provided  by 
the  DD-CM  algorithm,  we  would  not  have  been  able  to  isolate 
and  understand  these  different  sources  of  uncertainty  in  areas  2 
and  8. 


Next,  we  calculate  the  mean  magnetization  strength  for  all 
study  areas  using  the  best-fit  depths  and  magnetic  moments  re¬ 
covered  by  both  the  DD-CM  and  GD-VM  algorithms  (Table  1).  For 
this  calculation,  the  area  of  the  magnetic  anomalies  is  taken  as  the 
area  across  which  the  magnetic  field  B  >  Bmax/4}  where  Bmax  is  the 
peak  magnetic  field  in  a  given  study  area.  We  choose  to  approx¬ 
imate  the  thickness  of  the  magnetic  material  as  twice  the  aver¬ 
age  of  the  best-fit  magnetic  source  depths.  Most  values  cluster  in 
the  0.1 -0.3  A/m  range  found  by  Purucker  et  al.  (2012),  who  stud¬ 
ied  a  region  similar  to  area  1.  A  few  isolated  datasets  at  the  lowest 
altitudes  show  higher  magnetizations,  but  with  the  exception  of 
area  9,  average  magnetizations  for  a  particular  study  area  across  all 
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Fig.  14.  Area  7  best-fit  results  from  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  that  magnetometer  measurements  are  taken  from  K-MAG,  2009  day  123, 
altitude  31.7  km. 
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Fig.  15.  Area  8  best-fit  results  from  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  that  magnetometer  measurements  are  taken  from  K-MAG,  2009  day  123, 
altitude  31.4  km. 


datasets  are  close  to  the  0.1  -0.4  A/m  range  suggested  by  Wieczorek 
and  Weiss  (2010)  for  mare  basalts  and  mafic  impact  melt  breccias. 
Even  for  area  9,  only  the  GD-VM  average  magnetization  is  outside 
this  range. 

In  most  cases  the  GD-VM  algorithm  returns  higher  values  for 
magnetizations  than  the  DD-CM  algorithm.  This  is  a  result  of  the 
shallower  depths  returned  by  the  GD-VM  algorithm,  which  arise 
because  of  the  sheet-like  nature  of  the  magnetic  source,  compared 
to  the  dipolar  model  in  the  DD-CM  algorithm.  In  general  the  depth 
of  the  magnetic  source  bodies  is  not  well  constrained  from  our 
results,  due  to  the  possibility  of  non-uniquely  trading  depth  with 
moment. 


Finally,  best-fit  paleopoles  for  all  datasets  are  calculated  from 
the  best-fit  magnetization  directions  (Butler,  1998)  (Fig.  21, 
Table  3).  The  uncertainties  derived  from  the  methods  described  in 
the  Methods  section  yield  a  circular  error  ellipse  for  each  area  that 
becomes  distorted  when  its  paleopole  is  calculated.  Hence  we  cal¬ 
culate  the  two  different  semi-axes  of  the  ellipse  along  the  great 
circle  path  from  the  site  to  pole  (dp)  and  the  semi-axes  of  the  el¬ 
lipse  perpendicular  to  that  path  (dm)  (Butler,  1998).  For  area  2,  we 
also  show  the  paleopoles  from  the  two  local  error  minima  found 
in  the  DD-CM  method,  particularly  because  one  falls  very  close  to 
the  paleopole  for  area  1  (within  its  uncertainty),  while  the  other 
paleopoles  from  areas  3  and  6  (within  the  uncertainty  of  6). 
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Fig.  16.  Area  9  best-fit  results  from  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  that  magnetometer  measurements  are  taken  from  K-MAG,  2009  day  123, 
altitude  28.3  km. 
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Fig.  17.  Area  10  best-fit  results  from  both  inversion  algorithms.  Figure  details  are  as  in  Fig.  8,  except  that  magnetometer  measurements  are  taken  from  K-MAG,  2009  day 
123,  altitude  34.2  km. 
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Fig.  18.  RMS  error  map  (nT)  for  the  magnetization  direction  for  all  datasets  for  SPA  study  areas  1-5,  obtained  using  the  DD-CM  algorithm.  Southern  hemisphere  is  positive 
inclination;  the  whole-sphere  Mollweide  projection  is  centered  at  0°  inclination  and  0°  declination.  The  solid  white  contour  indicates  one  angular  standard  deviation  of 
dispersion  from  the  best-fit  solution  (red  dots,  Table  1),  from  Monte  Carlo  simulations  of  the  addition  of  1  nT  Gaussian  noise  to  observations,  using  the  DD-CM  algorithm. 
The  larger,  dashed  white  contour  indicates  the  minimum  error  solution  plus  1  nT,  a  measure  of  uncertainty  using  an  arbitrary  threshold  equal  to  the  measurement  uncertainty 
(not  used  in  our  final  analysis).  The  error  from  Monte  Carlo  simulations  of  the  effects  of  displacing  the  nominal  source  dipole  locations  by  0.5°  is  not  shown,  and  nor  is 
the  95%  confidence  interval  from  averaging  different  altitude  data  sets.  (For  interpretation  of  the  references  to  colour  in  this  figure  legend,  the  reader  is  referred  to  the  web 
version  of  this  article.) 


3.2.  Comparison  with  other  magnetic  paleopoles 

In  Fig.  22b-d  we  show  paleopoles  from  isolated  mag¬ 
netic  anomalies  studied  by  Takahashi  et  al.  (2014)  and 
Arkani-Hamed  and  Boutin  (2014),  as  well  as  paleopoles  inferred 
from  the  magnetization  of  the  crust  at  craters  (Arkani-Hamed 
and  Boutin,  2014).  The  magnetized  crust  paleopoles  from  Arkani- 
Hamed  and  Boutin  (2014)  are  the  means  of  the  values  listed  in 
their  Table  2.  None  of  the  anomalies  studied  by  these  groups  were 
inside  SPA.  To  compare  our  results  with  these  three  datasets,  we 
reverse  all  of  our  paleopoles  into  the  same  southern  hemisphere 
(Fig.  22a),  as  in  Takahashi  et  al.  (2014).  Here  we  have  dropped 
the  paleopole  found  from  the  mean  of  all  data  sets  at  area  2,  and 
only  show  the  paleopoles  from  the  two  local  error  minima  in  the 
DD-CM  algorithm. 

Takahashi  et  al.  (2014)  found  paleopoles  that  cluster  into  two 
groups,  one  near  the  present  pole,  and  another  at  mid-latitudes. 
We  find  a  wider  dispersion  in  the  distribution  of  our  paleopoles, 
but  we  do  find  some  clustering  at  the  present  pole  (areas  5,  7  and 
10),  as  well  as  some  near  their  mid-latitude  cluster  (areas  3,  6, 
and  8).  The  remaining  poles  are  not  easily  assigned  to  either  of 
these  clusters.  However,  overall,  the  paleopoles  we  find  do  seem 
to  avoid  longitudes  on  the  farside  (the  bottom  half  of  the  sphere  in 
Fig.  22a),  and  latitudes  <  30°.  None  of  the  paleopoles  from  Arkani- 
Hamed  and  Boutin  (2014)  show  obvious  correlation  with  any  of  the 
clustering  found  here  or  in  Takahashi  et  al.  (2014). 

Finally,  we  note  that  many  of  the  paleopoles  are  substantially 
separated  by  their  error  ellipses,  such  that  it  is  unlikely  that  only 


one  or  two  paleopoles  could  be  used  to  explain  all  of  the  anoma¬ 
lies,  without  severely  affecting  their  RMS  error  values.  This  can 
also  be  visualized  by  examining  the  error  spaces  in  Figs.  18  and 
19;  because  most  anomalies  are  at  similar  latitudes  and  longitudes 
(most  within  ~30°),  the  diversity  in  paleopole  locations  is  mostly 
determined  by  the  diversity  in  magnetization  directions.  Hence, 
the  diversity  in  RMS  error  minima  in  Figs.  18/19  also  graphically 
illustrates  the  diversity  in  paleopole  locations. 

In  the  next  section,  we  address  possible  origins  for  the  dis¬ 
persion  seen  in  our  study,  and  test  some  of  the  hypotheses  for 
anomaly  formation. 

4.  Discussion 

The  wide  variation  in  magnetic  paleopoles  derived  from  SPA’s 
magnetic  anomalies  presents  a  puzzle.  In  Sections  4.1 -4.4,  we  dis¬ 
cuss  the  implications  of  the  diverse  paleopole  locations  for  the 
formation  of  these  anomalies,  the  history  of  true  polar  wander,  and 
the  nature  of  the  lunar  dynamo. 

4.1  Magnetization  by  impact  processes 

A  long-standing  hypothesis  is  that  the  strongest  lunar  magnetic 
anomalies  are  genetically  related  to  the  antipodes  of  the  Imbrium, 
Orientale,  Serenitatis  and  Crisium  basins  (Hood  and  Williams, 
1989;  Hood  et  al.,  2001;  Lin  et  al.,  1998).  Compression  and  ampli¬ 
fication  of  the  IMF  by  impact-produced  plasma  may  be  strongest 
at  the  basin  antipode,  where  impact  ejecta  may  also  preferentially 
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Fig.  19.  RMS  error  map  (nT)  for  the  magnetization  direction  for  all  datasets  for  SPA  study  areas  6-10,  obtained  using  the  DD-CM  algorithm.  Description  is  as  in  Fig.  18. 


collect  (Hood  and  Artemieva,  2008).  The  antipodes  of  the  Imbrium 
and  Serenitatis  basins  are  close  to  our  study  areas  (Fig.  1),  suggest¬ 
ing  this  process  may  be  responsible  for  forming  the  anomalies  we 
have  examined  here.  If  true,  antipodal  ejecta  may  become  magne¬ 
tized  via  either  thermo-remanent  magnetization  (TRM)  or  shock- 
remanent  magnetization  (SRM).  Below  we  assess  these  two  possi¬ 
bilities. 

4  A  A.  TRM  in  an  impact-produced  field 

If  the  magnetization  were  produced  by  a  TRM,  hot  ejecta 
would  cool  in  the  presence  of  transient  IMF-amplified  fields  that 


would  last,  at  most,  for  one  day.  Assuming  a  thermal  diffusivity 
of  10_6m2/s,  the  thermal  cooling  timescale  for  one  day  is  <1  m. 
Therefore,  TRM  would  be  restricted  to  <1  m  of  material.  The  mag¬ 
netic  moments  we  find  across  all  SPA  study  areas  range  between 
1013Am2  and  1014Am2.  Using  source  body  horizontal  extents  from 
Fig.  1  (black  dashed  boxes),  depths  of  ~1  m  lead  to  high  TRM  mag¬ 
netizations  ranging  between  103  and  104A/m.  These  values  are  3- 
4  orders  of  magnitude  higher  than  samples  recovered  by  Apollo 
missions.  Additionally,  the  top  ~1  m  would  have  been  completely 
overturned  and  demagnetized  in  the  time  since  the  antipodal  im¬ 
pact  ~4  billion  years  ago  (Arnold,  1975).  Therefore,  TRM  from 
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impact-related  fields  is  implausible,  in  agreement  with  Hood  and 
Artemieva  (2008). 

4.12.  SRM  in  an  impact-produced  field 

Alternately,  impact  shock  pressures  from  the  deposition  of 
ejecta  at  the  antipode  may  create  SRM  in  the  ejecta  deposit.  How¬ 
ever,  to  allow  unidirectional  SRM  to  be  stably  imparted  to  a  large- 
scale  (e.g.  >4000 km2  for  areas  1  and  2;  >30,000 km2  for  areas 
5-8)  geologic  formation,  the  shock  waves  must  pass  through  the 
rock  when  the  entire  region  is  at  rest.  Shock  waves  begin  propa¬ 
gating  through  the  rock  instantly  upon  impact,  yet  at  that  moment 
the  ejecta  is  still  traveling  at  2.0-2.4km/s  (Hood  and  Artemieva, 
2008).  Therefore,  the  ejecta  at  the  antipode  will  change  orientation 
as  it  comes  to  rest,  after  any  SRM  has  formed,  such  that  the  total 
remanence  of  the  ejecta  deposit  will  be  randomized  and  nulled. 
Hood  and  Artemieva  (2008)  also  rule  out  SRM  in  the  ejecta,  using 
a  different  line  of  reasoning. 

It  is  possible  that  the  underlying  rock  could  be  shocked  by 
the  impacting  ejecta  and  acquire  a  unidirectional  SRM  (Hood  and 
Artemieva,  2008).  However,  we  suggest  that  there  are  two  prob¬ 
lems  with  this  hypothesis.  Firstly,  one  would  expect  nearly  verti¬ 
cal  local  magnetization  if  the  SRM  producing  field  was  compressed 
by  solar  wind  plasma  (see,  for  example,  field  lines  in  Figs.  8-9  of 
Hood  and  Artemieva,  2008).  Instead,  only  one  magnetic  anomaly, 
area  4,  has  a  nearly  vertical  magnetization,  and  the  rest  show  a 
preference  for  low  inclinations,  if  anything  (the  rest  are  at  least 
>  33°  away  from  the  vertical,  Fig.  20). 

Secondly,  deposition  of  impact  ejecta  would  take  place  over 
a  short  period  of  time,  over  which  the  ambient  field  direction 
is  likely  to  be  nearly  constant.  There  are  only  two  basin  an¬ 
tipodes  within  SPA,  which  implies  that  the  magnetic  anomalies 
should  have  one  of  two  magnetization  directions.  Roughly,  areas  1- 
3  (northwest  cluster),  4,  and  9  are  antipodal  to  Imbrium,  while  ar¬ 
eas  5-8  (eastern  stripes),  and  10  are  antipodal  to  Serenitatis.  First, 
we  consider  the  five  anomalies  at  the  Imbrium  antipode.  Here  we 
find  that  areas  1  and  3  have  magnetization  vectors  109°  apart 
(Fig.  21).  The  magnetization  at  area  4  is  95°  from  that  at  area  1  and 
70°  from  area  3.  The  magnetization  at  area  9  (Mare  Ingenii)  is  104° 
from  area  1,  140°  from  area  3,  and  128°  from  area  4.  All  of  these 
separations  are  well  outside  the  error  ellipses.  The  only  two  clus¬ 
ters  are  the  directions  of  areas  1  and  2,  and  possibly  area  3  with 
the  second  error  minimum  obtained  at  area  2.  In  sum,  there  are 
at  least  four  widely  separated  magnetization  directions  near  the 
Imbrium  antipode,  all  well  separated  by  their  error  ellipses.  It  is 
plausible  the  similar  directions  at  areas  1  and  2  suggest  this  pair 
was  magnetized  by  deposition  of  impact  ejecta,  but  it  would  not 
explain  the  magnetization  of  the  other  areas. 

At  the  five  anomalies  near  the  Serenitatis  antipode,  areas  7,  8, 
and  10  have  similar  directions,  but  5  and  6  are  well  separated  from 
each  other  and  from  7,  8,  and  10  (Fig.  21).  Therefore,  we  conclude 
again  that  it  is  not  likely  that  the  Serenitatis  impact  is  responsible 
for  magnetizing  all  of  these  anomalies. 

In  sum,  the  diversity  of  directions  argues  against  the  SRM  hy¬ 
pothesis,  or  at  least  allows  only  a  subset  of  the  geographically  clus¬ 
tered  anomalies  to  be  due  to  SRM.  Similar  arguments  can  be  ap¬ 
plied  to  ruling  out  SRM  from  surface  seismic  waves.  Considering 
all  of  the  observations  above,  we  conclude  that  at  least  some  of 
the  magnetic  source  bodies  in  SPA  were  magnetized  in  a  lunar  dy¬ 
namo  field,  rather  than  a  field  associated  with  impact  events. 

4.2.  Magnetized  South  Pole-Aitken  basin  ejecta  or  volcanic  bodies? 

Wieczorek  et  al.  (2012)  proposed  that  material  from  the  SPA 
impactor  might  be  the  source  of  many  of  the  SPA  basin’s  mag¬ 
netic  anomalies,  and  even  other  anomalies  across  the  Moon.  Un¬ 
der  this  hypothesis,  hot  iron-enriched  material  from  the  SPA  im¬ 


pactor  acquired  a  TRM  in  a  dynamo  field.  However,  if  SPA  ejecta 
landed  hot  and  cooled  in  a  dynamo  field,  the  resulting  anoma¬ 
lies  should  all  have  the  same  magnetization  direction  (but  see 
Section  4.3.4).  Instead,  the  diverse  magnetization  directions  (see 
Section  4.1.2,  above)  suggest  that  they  were  magnetized  at  differ¬ 
ent  epochs. 

The  last  remaining  viable  hypothesis  is  that  magnetic  anoma¬ 
lies  in  the  northern  SPA  basin  formed  as  a  result  of  magnetized 
sub-surface  dikes  (Purucker  et  al.,  2012).  Dikes  forming  over  long 
time  periods  would  permit  different  magnetization  directions  dur¬ 
ing  different  magnetic  epochs.  The  cluster  of  three  paleopoles  close 
to  the  present  pole  (areas  5,  7  and  10,  Fig.  22a)  would  suggest 
a  traditional  axial-aligned  dynamo  magnetized  these  dikes  when 
the  Moon  was  in  its  present  orientation.  However,  the  diversity 
in  paleopoles  seen  is  still  enigmatic  (see  Section  4.3).  Further,  if 
the  dike  hypothesis  is  true,  this  implies  that  the  dikes  near  SPA 
are  special  in  some  way,  since  the  nearside  of  the  Moon,  covered 
much  more  extensively  by  volcanism  (and  presumably  associated 
with  subsurface  dikes),  shows  no  magnetic  structures  like  those  at 
SPA.  Andrews-Hanna  et  al.  (2014)  reported  linear  gravity-gradient 
anomalies  that  may  be  dikes,  but  so  far,  no  obvious  correlation  be¬ 
tween  these  structures  and  magnetic  anomalies  have  been  found. 
We  note  that  Gong  and  Wieczorek  (2016)  find  a  correlation  be¬ 
tween  magnetization  and  gravity  anomalies  in  some  locations. 

4.3.  Possible  explanations  for  diverse  magnetic  field  directions 

The  above  arguments  suggest  the  magnetizing  fields  for  many 
of  the  SPA  anomalies  arose  from  a  dynamo.  However,  they  do  not 
offer  an  obvious  explanation  for  the  diversity  in  field  directions. 
Due  to  the  small  size  of  the  lunar  core,  and  the  rapid  decay  of 
magnetic  quadrupole  and  higher  terms  as  a  function  of  distance, 
it  is  likely  that  the  ancient  lunar  dynamo  was  dominantly  dipolar 
at  the  surface  (Weiss  and  Tikoo,  2014).  If  the  dipole  was  aligned 
with  the  Moon’s  spin  axis,  the  magnetization  directions  contain  in¬ 
formation  about  the  Moon’s  paleopole.  The  diverse  paleopoles  in 
Fig.  21  seem  to  imply  large  amounts  of  true  polar  wander 
(Goldreich  and  Toomre,  1969;  Runcorn,  1983).  While  the  diversity 
in  Fig.  21  is  surprising,  in  comparison  with  the  results  of  (Arkani- 
Hamed  and  Boutin,  2014;  Takahashi  et  al.,  2014)  (Fig.  22a),  the  di¬ 
versity  of  paleopoles  seen  in  our  findings  is  still  greater.  Below  we 
consider  some  possible  explanations  for  the  diverse  paleopole  lo¬ 
cations. 

4.3.1.  Non-axially  aligned  dipoles  and  impact-induced  dynamos 

Currently,  our  understanding  of  the  nature  of  the  lunar  dynamo 

is  limited.  Paleomagnetic  studies  favor  a  dynamo  that  existed  from 
approximately  4.2-3.3  Ga  ago.  Mechanisms  for  sustaining  such  a 
long-lived  core  dynamo  are  uncertain,  with  recent  proposals  for 
dynamos  driven  by  mechanical  stirring  from  impacts  and  preces¬ 
sion  (Dwyer  et  al.,  2011;  Le  Bars  et  al.,  2011;  Weiss  and  Tikoo, 
2014).  However,  it  is  not  known  if  these  dynamos  produce  the 
same  field  organization  as  the  Earth’s  dynamo.  It  may  be  possi¬ 
ble  that  these  exotic  dynamos  exhibit  more  variable  dipole  axis  di¬ 
rections,  which  could  explain  the  diversity  of  paleopoles  seen.  Un¬ 
fortunately,  more  work  modeling  small  dynamos  like  that  of  the 
Moon  is  needed  before  evaluating  this  hypothesis  further. 

4.3.2.  Secular  variation 

Another  possible  origin  for  the  diverse  paleopole  locations  is 
geomagnetic  secular  variation  of  a  dynamo  that  is  on  average 
aligned  with  the  lunar  spin  axis.  The  Earth’s  spin  axis  is  presently 
11°  away  from  its  magnetic  dipole  axis,  and  it  is  plausible  that  the 
Moon’s  ancient  dipole  axis  was  also  not  exactly  aligned  with  its 
spin  axis.  It  is  also  possible  that  secular  drift  might  be  larger  on  a 
body  with  a  small  core  such  as  the  Moon.  However,  some  of  the 
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Fig.  20.  Best-fit  magnetization  directions  for  study  areas  defined  in  Fig.  1,  for  ail  datasets  in  Figs.  2  and  3.  Directions  derived  from  different  altitudes  and  different  spacecraft 
at  the  same  study  area  are  represented  by  the  same  color.  Positive  inclinations  are  in  the  lower  hemisphere  (filled  circles).  (A  &  B)  Results  for  the  DD-CM  algorithm. 
Ellipses  in  A  represent  one  angular  standard  deviation  of  dispersion  from  1  nT  noise  simulations.  Ellipses  in  B  represent  one  angular  standard  deviation  of  dispersion  from 
simulations  of  source  dipoles  displaced  from  their  nominal  positions  (0.5°),  (C)  Results  from  the  GD-VM  algorithm,  including  results  from  regressions  to  merged  data  at  all 
altitudes  shown  in  Fig.  2  and  Fig.  3,  within  a  given  study  area  (stars). 


paleopole  locations  at  SPA  are  over  ~75°  apart  (after  accounting 
for  the  possibility  of  reversals,  Fig.  22a),  and  over  the  timescale  of 
magnetic  anomaly  formation,  the  mean  dipole  orientation  might 
average  to  be  aligned  with  the  spin  axis.  Presently,  we  have  little 
information  available  about  the  lunar  dynamo  to  further  evaluate 
the  role  of  secular  variation  in  explaining  the  diversity  of  paleopole 
directions. 


4.3.3.  True  polar  wander  over  long  timescales 

The  large  amount  of  polar  wander  implied  by  the  paleopoles 
is  difficult  to  reconcile  with  other  geophysical  constraints  for  the 
orientation  history  of  the  Moon.  Currently  there  are  two  compre¬ 
hensive  studies  that  derive  estimates  of  the  degree  of  polar  wan¬ 
der  on  the  Moon  due  to  long-term  changes  in  the  Moon’s  mo¬ 
ments  of  inertia  (Garrick-Bethell,  2016).  The  first  uses  the  shape 
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Fig.  21.  SPA  magnetic  paleopoles.  Paleopoles  are  the  mean  of  the  combined  results  from  the  DD-CM  and  GD-VM  algorithms,  from  all  directions  (altitudes)  shown  in  Fig.  20. 
Ellipses  represent  the  95%  confidence  interval  obtained  when  calculating  the  mean  direction  from  the  combined  DD-CM  and  GD-VM  results,  plus  the  dispersion  from  Monte 
Carlo  simulations  of  time-dependent  fields  and  dipoles  displaced  from  their  nominal  locations  (see  Methods).  Error  ellipses  are  not  calculated  for  the  two  different  minimum 
error  solutions  for  area  2  from  the  DD-CM  algorithm. 


and  gravity  of  the  Moon  exterior  to  large  basins  to  establish  the 
earliest  orientation  of  the  Moon  (Garrick-Bethell  et  al.,  2014).  The 
authors  find  the  present  lunar  pole  has  changed  by  ~36°  from 
its  earliest  axis,  but  it  does  not  coincide  with  any  of  the  pale¬ 
opole  clusters  we  find  (Fig.  22a).  The  second  study  uses  polar 
hydrogen  deposits  as  a  constraint  on  the  history  of  polar  wan¬ 
der,  and  infers  that  up  to  ~10°  of  paleolatitude  change  has  oc¬ 
curred  (Siegler  et  al.,  2016).  However,  they  cannot  readily  pro¬ 
duce  the  large  paleopole  changes  implied  by  the  anomalies  stud¬ 
ied  here.  In  summary,  the  limited  number  of  available  geophysical 
models  that  estimate  lunar  polar  wander  cannot  produce  the  di¬ 
versity  and  large  magnitude  of  paleopole  changes  required  to  ex¬ 
plain  our  observations.  Interestingly,  we  find  that  all  of  our  pale¬ 
opoles  are  >  35°  from  the  present  equator  (Fig.  22a,  accounting  for 
reversals),  which  may  be  due  to  the  difficulty  in  producing  the 
large  required  changes  in  the  Moon’s  moments  of  inertia. 


4.3.4.  True  polar  wander  due  to  SPA  formation 

Very  large  initial  changes  in  the  Moon’s  moments  of  inertia  due 
to  SPA’s  crater  might  have  produced  large  amounts  of  polar  wan¬ 
der.  Eventually,  these  changes  must  have  subsided  over  millions 
to  billions  of  years,  since  the  gravity  signature  observed  today  is 
muted  (Zuber  et  al.,  2013).  If  hot  material  was  cooling  in  the  pres¬ 
ence  of  a  dynamo  field  throughout  the  Moon’s  reorientation,  its 
magnetization  could  in  principle  capture  multiple  lunar  orienta¬ 
tions,  thereby  producing  the  diversity  of  paleopole  locations.  Some 
of  this  hot  material  could  be  iron-rich  material  from  the  SPA  im- 
pactor  (Wieczorek  et  al.,  2012).  Following  impact,  the  free  preces¬ 
sion  damping  time  of  the  Moon  is  ~2  x  105  years  (Peale,  1976; 
Williams  et  al.,  2001),  assuming  a  semimajor  axis  of  30  Earth  radii, 
dissipation  quality  factor  Q.=  50  and  degree-2  love  number  k2  =  0.1 
(which  we  assume  to  be  representative  of  the  Moon  when  SPA 
formed).  The  length  scale  for  cooling  in  ~2  x  105  years  is  ~3  km 
(assuming  a  thermal  diffusivity  of  10-6  m2  s_1),  which  is  small 
compared  to  the  scale  of  the  anomalies  we  observe,  but  perhaps 
not  so  much  as  to  preclude  recording  a  measurable  magnetization. 


However,  using  a  simple  model  for  reorientation,  we  find  that 
density  anomalies  produced  by  SPA’s  crater  do  not  produce  pa¬ 
leopoles  that  overlap  with  those  of  its  magnetic  anomalies.  In 
our  model,  we  replace  SPA’s  gravity  potential  inside  the  outer  to¬ 
pography  rim  with  values  between  -3  times  the  maximum,  and 
+  2  times  the  maximum  of  the  present  day  potential,  in  incre¬ 
ments  of  0.25  times  the  maximum  potential.  The  inertia  tensor 
and  paleopoles  were  then  calculated  for  each  of  these  cases  us¬ 
ing  the  resulting  globally  calculated  degree-2  spherical  harmonics 
(cf.  Garrick-Bethell  et  al.,  2014).  These  positive  and  negative  gravity 
values  represent  very  different  models  of  SPA’s  effects,  but  illus¬ 
trate  the  range  of  paleopoles  that  are  possible.  We  find  the  max¬ 
imum  extent  of  true  polar  wander  from  SPA’s  formation  passes 
close  to  the  paleopoles  from  areas  7  and  10  (Fig.  22a,  the  negative 
anomaly  path  approaches  the  center  of  SPA  as  the  anomaly  size 
grows,  as  expected).  However,  these  anomalies  are  already  close  to 
the  present  pole,  and  do  not  need  to  be  explained  by  polar  wander. 
We  do  find  that  the  limits  of  polar  wander  approach  areas  4  and  8, 
such  that  this  process  could  plausibly  explain  some  of  these  mag¬ 
netic  paleopoles.  Unfortunately,  the  paleopoles  from  areas  1,  2,  3,  5, 
and  9,  and  to  a  lesser  extent  area  6,  do  not  come  close  to  the  pale¬ 
opoles  produced  by  SPA.  Thus,  cooling  of  SPA  ejecta  deposits  over 
short  timescales  (~105  years),  or  even  dikes  that  formed  when  the 
Moon  resided  at  any  of  the  SPA-produced  paleopoles,  cannot  fully 
explain  the  observed  diversity  of  magnetic  paleopoles  we  find. 

There  are  many  unknowns  in  modeling  the  paleopoles  al¬ 
lowed  by  SPA’s  formation.  For  example,  reorientation  and  mag¬ 
netic  anomaly  formation  depends  on  the  SPA  impact’s  effect  on  dy¬ 
namo  operation  (Arkani-Hamed  and  Olson,  2010a,  2010b),  the  in¬ 
ertia  tensor  of  the  Moon  just  before  SPA  formation  and  the  length 
scale  (and  thereby  cooling  timescale)  of  the  materials  making  up 
the  anomalies.  Furthermore,  the  spin  vector  of  the  Moon  will  be 
freely  precessing  around  its  angular  momentum  vector  during  the 
free  precession  damping  timescale.  If  the  dynamo  after  SPA  impact 
retained  its  alignment  along  the  Moon’s  angular  momentum  vec¬ 
tor,  this  precession  could  broaden  the  range  of  paleopoles  permit¬ 
ted  (essentially  accessing  a  range  of  paleopoles  around  the  pale- 
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Fig.  22.  Paleopoles  from  SPA  anomalies  compared  with  published  paleopoles  from  other  anomalies.  All  figures  show  the  southern  hemisphere  (all  points  are  southern 
latitudes).  (A)  Our  results  from  Fig.  21,  with  all  north  poles  reversed  into  the  southern  hemisphere.  The  paleopole  obtained  from  the  mean  of  all  data  at  area  2  (see  Fig. 
21)  has  been  omitted  in  this  figure.  Red  and  blue  triangles  are  possible  lunar  paleopoles  resulting  from  SPA’s  gravity  anomaly  immediately  after  its  formation,  assuming 
an  extreme  range  of  density  anomaly  models;  red  (blue)  points  are  positive  (negative)  density  anomalies  (see  Section  4.3.4).  SPA’s  center  (a  limiting  paleopole  for  very 
negative  density  models  of  SPA)  is  shown  at  (-53.2°  S,  191°  E)  (black  triangle).  The  magenta  triangle  shows  the  paleopole  dervied  from  the  tidal  component  of  the  Moon’s 
topography,  outside  of  large  basins  (Garrick-Bethell  et  al.,  2014).  (B)  Paleopoles  from  the  11  anomalies  reported  by  Takahashi  et  al.  (2014).  No  anomalies  are  within  SPA.  Blue 
(red)  points  represent  inversions  from  Kaguya  (Lunar  Prospector)  data.  In  many  cases,  multiple  points  represent  paleopoles  from  single  sites.  Stars  represent  the  means  of 
the  two  clusters.  (C)  Paleopoles  from  10  isolated  magnetic  anomalies  reported  by  Arkani-Hamed  et  al.  (2014).  (D)  Paleopoles  from  crustal  magnetization  at  20  impact  craters 
reported  by  Arkani-Hamed  et  al.  (2014).  (For  interpretation  of  the  references  to  color  in  this  figure  legend,  the  reader  is  referred  to  the  web  version  of  this  article.) 
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opoles  shown  in  Fig.  22a,  with  the  range  depending  on  the  preces¬ 
sion  angle).  Therefore,  we  cannot  definitively  rule  out  polar  wander 
processes  as  the  origin  for  some  of  the  diverse  paleopole  locations 
at  SPA. 

A  variant  of  the  hypothesis  of  iron-rich  SPA  ejecta  (Wieczorek 
et  al.,  2012)  is  the  cooling  of  iron-rich  “sesquinary”  (Zahnle  et  al., 
2008)  impactors  formed  by  the  SPA  impact.  Ejected  into  orbit  im¬ 
mediately  after  impact,  studies  on  various  planetary  bodies  show 
that  these  impactors  can  return  to  a  body  from  104  years  (Nayak 
et  al.,  2016)  to  106  years  (Gladman  et  al.,  1995)  post-impact,  at 
approximately  escape  velocity  (Nayak  and  Asphaug,  2016).  With 
reimpacts  spread  across  this  timescale,  reaccreting  iron-rich  mate¬ 
rial  originally  from  SPA,  either  still  hot  from  the  impact  or  heating 
upon  re-impact,  could  record  a  diverse  set  of  orientations  as  the 
Moon  reorients  in  response  to  moment  of  inertia  changes.  Across 
this  timescale,  impact  locations  become  random;  such  sesquinary 
magnetism  caused  by  these  iron-rich  impactors  would  be  widely 
distributed  around  the  Moon. 

5.  Conclusions 

Using  two  different  inversion  methods,  we  find  diverse  direc¬ 
tions  of  magnetization  among  magnetic  anomalies  in  the  north¬ 
ern  SPA  basin.  The  diverse  directions  help  rule  out  impact-related 
fields  as  their  only  origin.  Intrusive  bodies,  such  as  the  dikes  pro¬ 
posed  by  (Purucker  et  al.,  2012),  are  a  plausible  explanation.  The 
diverse  paleopole  locations  could  imply  large  amounts  of  true  po¬ 
lar  wander,  but  true  polar  wander  inferred  independently  from 
gravity  (Garrick-Bethell  et  al.,  2014)  and  hydrogen  deposits  (Siegler 
et  al.,  2016)  implies  more  modest  changes  in  the  Moon’s  orienta¬ 
tion.  The  diverse  directions  argue  against  the  hypothesis  that  they 
all  formed  from  iron-rich  SPA  ejecta  that  cooled  in  a  dynamo  field 
(Wieczorek  et  al.,  2012).  A  simple  gravity  anomaly  model  for  large 
amounts  of  SPA-produced  reorientation  fails  to  explain  at  least  five 
of  the  paleopoles  observed,  but  many  unknowns  remain  in  mod¬ 
eling  this  process.  Some  SPA  ejecta  may  have  produced  iron-rich 


sesquinary  impactors  that  landed  across  the  Moon  and  recorded 
orientation  changes  as  they  cooled  in  a  dynamo  field.  A  dynamo 
that  was  not  aligned  with  the  lunar  spin  axis  remains  a  plausi¬ 
ble  explanation  for  all  observations,  but  gaps  remain  in  our  ability 
confirm  this  hypothesis.  The  wide  variety  of  viable  hypotheses  and 
large  number  of  unknowns  highlight  the  complexity  of  interpreting 
the  origins  of  lunar  magnetic  anomalies  and  their  paleopoles. 
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Appendix  A.  Comparison  of  magnetic  field  observations  with 
best-fit  models  for  all  datasets 

We  use  45  datasets  for  the  analysis  in  this  work,  the  total  field 
plots  of  which  are  shown  in  Fig.  2  (areas  1-5)  and  Fig.  3  (areas 
6-10).  In  Figs.  8-17  we  presented  example  magnetic  field  maps  for 
the  best-fit  source  magnetizations,  from  both  inversion  algorithms, 
for  areas  1-10,  respectively.  Here  we  present  the  magnetic  field 
maps  for  the  best-fit  source  magnetizations  for  the  35  remaining 
datasets. 
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Fig.  Al.  Area  1  best-fit  results  for  both  inversion  algorithms.  (Left  column)  K-MAG  observations  in  the  lunar  wake  on  2009  day  123,  compared  to  (Middle  column)  a  model 
with  manually  placed  dipoles  of  equal  magnetization  (DD-CM)  and  (Right  column)  a  grid-based  model  with  variable  magnetizations  (GD-VM).  Source  dipoles  in  the  DD-CM 
algorithm  are  approximately  located  at  maxima  of  the  observed  total  field  (black  dots).  White  dots  in  the  east  panel  are  locations  of  spacecraft  magnetometer  measurements. 
The  mean  measurement  altitude  is  31.6  km.  Longitude  (degrees  east)  is  indicated  on  the  x-axis. 
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Fig.  A2.  Area  1  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  31.1  km  and  the  observations  are  from 
1999  day  61. 
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Fig.  A3.  Area  1  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  21.4  km  and  the  observations  are  from 
1999  day  142. 
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Fig.  A4.  Area  2  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  43.5  km  and  the  observations  are  from 
2009  day  69. 
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Fig.  A5.  Area  2  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  38.5  km  and  the  observations  are  from 
2009  day  96. 
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Fig.  A6.  Area  2  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  31.1  km  and  the  observations  are  from 
1999  day  33. 
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Fig.  A7.  Area  3  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  42.8  km  and  the  observations  are  from 
2009  day  69. 
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Fig.  A8.  Area  3  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  37.6  km  and  the  observations  are  from 
2009  day  96. 
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Fig.  A9.  Area  3  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  29.6  km  and  the  observations  are  from 
1999  day  61. 
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Fig.  A10.  Area  4  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  12.1  km  and  the  observations  are  from 
2009  day  151. 
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Fig.  All.  Area  4  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  30.9  km  and  the  observations  are  from 
1999  day  115. 
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Fig.  A12.  Area  4  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  21.8  km  and  the  observations  are  from 
1999  day  142. 
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Fig.  A13.  Area  5  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  39.1  km  and  the  observations  are  from 
2009  day  96. 


Fig.  A14.  Area  5  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  13.2  km  and  the  observations  are  from 
2009  day  151. 
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Fig.  A15.  Area  5  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  22.3  km  and  the  observations  are  from 
1999  day  142. 
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Fig.  A16.  Area  6  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  39.1  km  and  the  observations  are  from 
2009  day  96. 
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Fig.  A17.  Area  6  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  32.8  km  and  the  observations  are  from 
1999  day  33. 
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Fig.  A18.  Area  6  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  32.5  km  and  the  observations  are  from 
1999  day  115. 
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Fig.  A19.  Area  6  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  22.6  km  and  the  observations  are  from 
1999  day  142. 


Fig.  A20.  Area  7  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  39.0  km  and  the  observations  are  from 
2009  day  96. 
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Fig.  A21.  Area  7  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  32.6  km  and  the  observations  are  from 
1999  day  33. 


Fig.  A22.  Area  7  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  32.4  km  and  the  observations  are  from 
1999  day  115. 
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Fig.  A23.  Area  8  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  39.4  km  and  the  observations  are  from 
2009  day  96. 
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Fig.  A24.  Area  8  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  13.0  km  and  the  observations  are  from 
2009  day  151. 
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Fig.  A25.  Area  8  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  32.7  km  and  the  observations  are  from 
1999  day  61. 


Fig.  A26.  Area  8  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  23.6  km  and  the  observations  are  from 
1999  day  142. 
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Fig.  A27.  Area  9  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al ,  except  that  the  mean  measurement  altitude  is  37.0  km  and  the  observations  are  from 
2009  day  96. 
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Fig.  A28.  Area  9  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  28.7  km  and  the  observations  are  from 
1999  day  61. 
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Fig.  A29.  Area  9  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  28.5  km  and  the  observations  are  from 
1999  day  88. 
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Fig.  A30.  Area  9  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  28.7  km  and  the  observations  are  from 
1999  day  115. 
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Fig.  A31.  Area  9  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  19.6  km  and  the  observations  are  from 
1999  day  142. 
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Fig.  A32.  Area  10  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  13.9  km  and  the  observations  are 
from  2009  day  151. 
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Fig.  A33.  Area  10  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  32.3  km  and  the  observations  are 
from  1999  day  61. 
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Fig.  A34.  Area  10  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  32.0  km  and  the  observations  are 
from  1999  day  88. 
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Fig.  A35.  Area  10  best-fit  results  for  both  inversion  algorithms.  Details  are  as  in  Fig.  Al,  except  that  the  mean  measurement  altitude  is  21.9  km  and  the  observations  are 
from  1999  day  142. 
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